Nuclear Force from Lattice QCD 
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Nucleon-nucleon (NN) potential is studied by lattice QCD simulations in the quenched approx- 
imation, using the plaquette gauge action and the Wilson quark action on a 32 4 (~ (4.4 fm) 4 ) 
lattice. A NN potential VnnM is denned from the equal-time Bethe-Salpeter amplitude with a 
local interpolating operator for the nucleon. By studying the NN interaction in the 1 So and 3 Si 
channels, we show that the central part of VW(r) has a strong repulsive core of a few hundred MeV 
at short distances (r < 0.5 fm) surrounded by an attractive well at medium and long distances. 
These features are consistent with the known phenomenological features of the nuclear force. 

PACS numbers: 12.38.Gc, 13.75.Cs, 21.30-Cb 



More than 70 years ago, Yukawa introduced the pion 
to account for the strong interaction between the nucle- 
ons (the nuclear force) Since then, enormous efforts 
have been devoted to understand the nucleon-nucleon 
(NN) interaction at low energies both from theoretical 
and experimental points of view. As shown in Fig[TJ 
phenomenological NN potentials are thought to be char- 
acterized by three distinct regions [3, Q : The long range 
part (r > 2 fm) is well understood and is dominated 
by the one pion exchange. The medium range part 
(1 fm < r < 2 fm) receives significant contributions from 
the exchange of multi-pions and heavy mesons (p, uj, and 
a). The short range part (r < 1 fm) is empirically known 
to have strong repulsive core Q, which is essential not 
only for describing the NN scattering data, but also 
for the stability and saturation of atomic nuclei, for de- 
termining the maximum mass of neutron stars, and for 
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FIG. 1: Three examples of the modern NN potential in the 
1 So (spin singlet and s-wave) channel: CD-Bonn 4], Reid93 
H and AV18 from the top at r = 0.8 fm. 



igniting the Type II supernova explosions 0. Although 
the origin of the repulsive core must be closely related to 
the quark-gluon structure of the nucleon, it has been a 
long-standing open questions in QCD [3]. 

In this Letter, we report our first serious attempt to 
attack the problem of nuclear force from lattice QCD 
simulations [13] • The essential idea is to define a NN 
potential from the equal-time Bethe-Salpeter (BS) ampli- 
tude of the two local interpolating operators separated by 
distance r [llj]. This type of BS amplitude has been em- 
ployed by CP-PACS collaboration to study the irn scat- 
tering on the lattice [13]. As we shall see below, our 
NN potential shows a strong repulsive core of about a 
few hundred MeV at short distances surrounded by an 
attraction at medium and long distances in the s-wave 
channel. 

Let us start with an effective Schrodinger equation ob- 
tained from the BS amplitude for two nucleons at low 
energies 12. 13j|: 
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where fi = ton/ 2 and E is the reduced mass of the nu- 
cleon and the non-relativistic energy, respectively. For 
the NN scattering at low energies, the non-local poten- 
tial U is represented as U(f, f) = Vnn<7, V)<5(r— f) with 
the derivative expansion [2j: Vnn = Vc(r) + Vr(r)Si2 + 
V LS (r)L ■ S + 0{V 2 ). Here S a = 3(3*1 • f){S 2 ■ f) - B x ■ a 2 
is the tensor operator with f = \r\/r, S the total spin 
operator, and L = —ir x V the relative angular mo- 
mentum operator. The central NN potential Vc(r), the 
tensor potential Vt (r) and the spin-orbit potential Vls ( r ) 
can be further decomposed into various spin-isospin chan- 
nels, e.g. V c (r) = Vg(r) + V§{r)a x ■ a 2 + VJ(r)n • f 2 + 
Vc T ( r ){&i ' ( ?2)(t ? i • In the phenomenological analy- 
sis of the NN scattering phase shift [3(] , the Schrodinger 
equation with a certain parametrization of Vnn is solved 
and compared with the data. On the other hand, if we 
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can calculate <j>(f) directly from lattice simulations for 
various E, Eq.|T]) can be used to define the non-local 
potential <7(r, f*) directly without recourse to the exper- 
imental inputs except for quark masses and the QCD 
scale parameter. In this Letter, instead of finding U by 
varying E, we take only the leading term in the deriva- 
tive expansion at low energies and extract the central 
potential Vfc(r) at fixed E through 



V c (r) = E- 



1 V 2 r/>(r) 
<f>(r) 



(2) 



On the lattice, <f>{r) with zero angular momentum (£ = 
0) is defined from the equal-time BS amplitude as 
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where we choose the local interpolating operator for the 
nucleon: N l a = e a bc ( t q a C^5T2q b ) q]f with a, b and c the 
color indices, a and (3 the Dirac indices, i and j the 
isospin indices, and C = 7472 the charge conjugation. 
r describes the spatial separation between the nucleons. 
Since we consider the NN scattering at low energies, we 
take only the upper components of N l a . The summation 
over the vector x projects out the state with zero total- 
momentum. The summation over discrete rotation TZ of 
the cubic group O projects out the representation 
which contains I — state and t > 4 states. The former 
can be singled out by selecting the lowest energy state 
with the procedure given in Eq.(j4])- The spin (isospin) 
projection is carried out by the operator P a (P T ); for 
example, P£p = (<T2)a^(= 5 a @) in the spin-singlet (spin- 
triplet) channel. The renormalization factor Z, which 
relates the BS amplitude on the lattice and that in the 
continuum, cancels out in Vc(r). 

The (f>(r — y — x) in Eq.([3]) is nothing but the prob- 
ability amplitude to find "nucleonlike" three quarks lo- 
cated at point x and another "nucleonlike" three quarks 
located at point y. In terms of the physical states, <p(f) 
contains not only the elastic amplitude NN— »NN but also 
the inelastic amplitudes such as NN — > 7rNN. However, 
at low energies below threshold, the inelastic part is spa- 
tially localized and does not affect the asymptotic form 
of 4>(f). Note also that a different choice of the nucleon 
interpolating operator modifies the relative weight of the 
elastic and inelastic amplitudes and thus leads to a dif- 
ferent NN potential. Nevertheless, they give the same 
scattering phase shift since the asymptotic form of </>(f) is 
independent of the interpolating operators. It is an open 
question at the moment how the short distant structure 
of Vnn to be shown below is affected by the change of the 
interpolating operator. For more details on these points, 
see 1131. 



In the actual simulations, Eq. ([3]) is obtained through 
the four point nucleon correlator, 



Fj<ii<i(x,y,t;to) = 



N l a {x,t)N J {y,t)J NN {t o ) 



K(x)N J (y) 



-E n (t-t ) 



(4) 



Here Jnn^o) is a source term located at t = to, which 
produces the nucleons with zero total momentum. To en- 
hance the ground state contribution of the NN system, we 
adopt the wall source, Jnn(*o) = P[ 3 P^JV^t Q )Af^(t ), 
where AT is obtained from N by replacing q by Q(to) = 
J2s 9(^> *o)- E n is the energy of the two-nucleon state \n) 
and A n (to) = (n\ J r NN(*o)|0}- Because of the finite lattice 
volume L 3 , the energy E takes discrete value and has a 
finite shift from the non-interacting case AE = 0(1/ X 3 ) 
to be determined from the simulations 

In this Letter, we focus on the spin-singlet and spin- 
triplet channels with zero orbital angular momentum. In 
the standard notation, the former (latter) corresponds to 
the 2s+1 £j =1 So (= 3 Si) channel, where s, £ and J denote 
the total spin, orbital angular momentum, and the total 
angular momentum of the two nucleons. The is the 
simplest channel where only the central potential Vc(r) 
contributes. On the other hand, there arises a mixing 
between the 3 Si and 3 Di channels because of the tensor 
force Vt (r) . In this case, one may define an effective cen- 
tral potential V^ ff (r) which consists of the bare central 
potential and the induced central potential by the 3 Di 
admixture 0- The definition in Eq.([2]) with <f>(r) being 
projected onto 1 So ( 3 Si) corresponds to the central po- 
tential (the effective central potential). 

To calculate 4>(r), we have carried out simulations on 
a 32 4 lattice in the quenched approximation. We em- 
ploy the plaquette gauge action with the gauge coupling 
(3 = 5.7 and the Wilson quark action. The lattice spacing 
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FIG. 2: The lattice QCD result of the radial dependence of 
the NN wave function at t— to — 6 in the 1 So and 3 Si channels. 
Inset shows the two-dimensional view in the x — y plane. 
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FIG. 3: The lattice QCD result of the central (effective cen- 
tral) part of the NN potential Vc(r) {V§ s {r)) in the ( 3 Si) 
channel for m^/mp = 0.595. The inset shows its enlargement. 
The solid lines correspond to the one-pion exchange potential 
(OPEP) given in Eq.©. 



determined from the p meson mass in the chiral limit is 
-i = 1.44(2) GeV (a ~ 0.137 fm) [15], which leads to 
the lattice size (4.4 fm) 4 . The hopping parameter is cho- 
sen to be k = 0.1665, which corresponds to to^ ~ 0.53 
GeV, to p ~ 0.89 GeV and m N ~ 1.34 GeV. We use the 
global heat-bath algorithm with overrelaxations to gener- 
ate the gauge configurations. After skipping 3000 sweeps 
for thermalization, 500 gauge configurations are collected 
with the interval of 200 sweeps. Results for lighter and 
heavier quark masses with higher statistics will be re- 
ported in [13 |. The Dirichlet (periodic) boundary con- 
dition for quarks is imposed in the temporal (spatial) 
direction. To avoid the boundary effect, the wall source 
is placed at t = to = 5 at which the Coulomb gauge fixing 
is made. The ground state saturation for t — to > 6 is 
checked by the effective mass of the two-nucleon system. 

Fig. [2] shows the lattice QCD result of the wave func- 
tion at the time-slice t — to = 6. They are normalized 
at the spatial boundary r = (32/2 = 16,0,0). All the 
data including the off- axis ones are plotted for r < 0.7 
fm, beyond which we plot only the data locating on the 
coordinate axes and their nearest neighbors. As is clear 
from Fig. [2] the wave function is suppressed at short dis- 
tance and have a slight enhancement at medium distance, 
which suggests that the NN system has a repulsion (at- 
traction) at short (medium) distance. 

Fig. [3] shows the central (effective central) NN poten- 
tial in the 1 So ( 3 Si) channel at t — 1 = 6. As for V 2 
in Eq. ([2]), we take the discrete form of the Laplacian 
with the nearest-neighbor points. E is obtained from 
the Green's function G(r; E) which is a solution of the 
Helmholtz equation on the lattice [l2|. By fitting the 
wave function c/>(r) at the points r = (10 — 16,0,0) and 
(10 - 16, 1, 0) by G(r; E), we obtain £(%) = -0.49(15) 



FIG. 4: t — to dependence of Vc(r) in the 1 So channel for 
several different values of the distance r. 



MeV and E( 3 Si) = -0.67(18) MeV. Namely, there is 
a slight attraction between the two nucleons in a finite 
box. To make an independent check of the ground state 
saturation, we plot the t dependence of Vc(r) in the 
channel at several distances r = 0,0.14,0.19,0.69,1.37 
and 2.19 fm in Fig. 0] The saturation indeed holds for 
t — to > 6 within errors. 

As anticipated from Fig. [21 Vc{r) and Vq S (t) have 
repulsive core at r < 0.5 fm with the height of about 
a few hundred MeV. Also, they have an attraction of 
about -(20-30) MeV at the distance 0.5 < r < 1.0 
fm. The solid lines in Fig. [3] show the one pion exchange 
contribution to the central potential calculated from 



VSir) 



gl N (n -t 2 )(cti - era) 
4?r 3 



(5) 



where we have used m, ~ 0.53 GeV and ton —1.34 GeV 
to be consistent with our data, while the physical value 
of the 7rN coupling constant is used, <7 2 .jv/(4 7r ) — 14.0. 
Even in the quenched approximation, the one pion ex- 
change is possible as the connected quark exchange be- 
tween the two nucleons. In addition, there is in principle 
a quenched artifact to the NN potential from the flavor- 
singlet hairpin diagram (the ghost exchange) between the 
nucleons [l6(. Its contribution to the central potential 

reads 0: V^r) = #^ (^f (l - ^) . 
Here g V N and too are the r]N coupling constant and a 
mass parameter of the ghost, respectively. The ghost 
potential has an exponential tail which dominates over 
the Yukawa potential at large distances. Its significance 
can be estimated by comparing the sign and the mag- 
nitude of e m7,r Vc(r) and e m " rr V^ ff (r) at large distances, 
because Vq (r) has an opposite sign between and 3 Si. 
Our present data shows no evidence of the ghost at large 
distances within errors, which may indicate g^jq <C g^N- 
Several comments are in order here: 
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1. The asymptotic wave function at low energy (E — > 0) 
is approximated as </> asy (r) = gM*i£+MM] _> r +ao ; where 
5o(fc) (°o) is the s-wave scattering phase shift (scattering 
length). From the zero of r4>asy(r), we find ao( 1 So) = 
0.066(22) fm and a ( 3 Si) = 0.089(27) fin under the as- 
sumption that E r~j —0.5 MeV is small enough [13]. The 
reason of having such a small ao is easily understood 
from the well-known formula in the Born approximation: 
ao — — ton Jq 00 Vc(r )r 2 dr, where (i) the volume factor 
r 2 dr hides the repulsive core at short distances, and (ii) 
ao is a subtle quantity which suffers from a large cancel- 
lation between repulsion and attraction. 

2. The above points (i) and (ii) also provide a reason why 
it is dangerous at the moment to compare ao obtained in 
unphysical quark masses with ao in experiments (l7j . A 
slight change of the depth and height of the potential by 
the change of the quark mass may affect the scattering 
length substantially. 

3. If the attraction becomes large enough at small m v , 
the system may have a bound state and ao changes sign 
from positive to negative. In the real world, this happens 
in 3 Si (the deuteron channel). Indeed, the inset of Fig. [3] 
shows that 3 Si has stronger attraction than 1 S - Note 
that this difference cannot be attributed to the ghost 
contribution, because it is repulsive (attractive) in 3 Si 
(%)■ 

4. Our preliminary result with the lighter quark mass 
at k = 0.1678 (m p ~ 0.36 GeV, m p ~ 0.84 GeV and 
tojv — 1.2 GeV) [lH shows that the height of the re- 
pulsive core increases by 40 % at the origin, while the 
minimum of the potential is shifted to r ~ 0.8 fm with 
approximately the same depth. 

In summary, we have studied the NN interaction using 
the lattice QCD simulations in the quenched approxima- 
tion on a (4.4 fm) 4 lattice with the quark mass corre- 
sponding to m w /m p — 0.595. We define a NN poten- 
tial with the use of the equal-time Bethe-Salpeter am- 
plitude for two local interpolating operators of the nu- 
cleon. The central (effective central) part of the poten- 
tial in the 1 So ( 3 Si) channel at low energies turns out 
to have a repulsive core at short distances surrounded 
by attractive well at medium and long distances. These 
properties are known to be the important features of the 
phenomenological NN potential. The long range tails of 
Vc(r) in both channels are consistent with the one-pion 
exchange within statistical errors. We did not find nu- 
merical evidence of the ghost exchange (quenched arti- 
fact) for m^/mp — 0.595. Preliminary calculation with 
the lighter quark mass (m^/mp = 0.433) shows that the 
short range repulsion becomes stronger and the range of 
the potential becomes wider. 

It would be quite interesting to derive the tensor and 
spin-orbit forces by making appropriate projections of 
the wave function. Studies of the hyperon-nucleon and 
hyperon-hyperon potentials, whose experimental infor- 
mation is currently very limited, are now under investi- 



gation: They are particularly important for the physics of 
hyper nuclei and neutron star core. To unravel the phys- 
ical origin of the repulsive core, we need further studies 
on its quark mass dependence, channel dependence and 
interpolating-operator dependence. Also, the full QCD 
simulations at the physical quark mass are necessary to 
be done in the future. 
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